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ABSTRACT 

In this letter we study the eccentricity evolution of a massive black hole (MBH) binary 
(MBHB) embedded in a rotating stellar cusp. Following the observation that stars on 
counter-rotating (with respect to the MBHB) orbits extract angular momentum from 
the binary more efficiently then their co-rotating counterparts, the eccentricity evo- 
lution of the MBHB must depend on the degree of co-rotation (counter-rotation) of 
the surrounding stellar distribution. Using an hybrid scheme that couples numerical 
three-body scatterings to an analytical formalism for the cusp-binary interaction, we 
verify this hypothesis by evolving the MBHB in spherically symmetric cusps with dif- 
ferent fractions of co-rotating stars. Consistently with previous works, binaries in 
isotropic cusps {J- = 0.5) tend to increase their eccentricity, and when T approaches 
zero (counter-rotating cusps) the eccentricity rapidly increases to almost unity. Con- 
versely, binaries in cusps with a significant degree of co-rotation {T > 0.7) tend to 
become less and less eccentric, circularising quite quickly for T approaching unity. 
Direct A''-body integrations performed to test the theory, corroborate the results of 
the hybrid scheme, at least at a qualitative level. We discuss quantitative differences, 
ascribing their origin to the oversimplified nature of the hybrid approach. 
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1 INTRODUCTION 

Since their theoretical pred ic tion in the early 80's 
iBegelman. Blandford fc Rees I Il980l ). massive black 
hole (MBH) binaries (MBHBs) forming in galactic nuclei 
following galaxy mergers, have been the main focus of 
many dynamical studies. In the con text of the hierarchica l 
formation of galactic structures (jWhite fc Rees I Il978l ) 
and the MBHs residing at their center, bound MBH 
binaries forming at parsec scales have to get rid of their 
orbital energy to reach the point where gravitational 
waves (GW) emission becomes efficient enough to drive 
their final coalescence . This is known as the 'last parsec 
problem' (jMilosavlievic fc MerrittlbOOll ). Interactions with 
ambient stars, abundant in dense nuclei, might provide 
the physical source of energy extraction by means of the 
slingshot mechanism. That is, in a strong three-body 
encounter, the light intruder star is ejected at infinity 
carrying away part of the orbital energy and angular 
mom entum of the massive binary jMikkola fc Valtonen I 
Il992l ). In the la st two dec a des se ve ral analytical and numer 
ical works (iQuinlan I Il996l: iMilosavlievic fc Merritt 
I2OOII : Hemsendorf. Sigurdsson fc Spurzem 
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IMilosavlievic fc Merritt 



|2004 



20031: 



iMakino fc Funato 



iBerczik et al 



Baumgardt, Gualandris fc P ortcgics Zwart 
I2OO6I: "iMerritt fc Szell 
iMerritt. Mikkola fc Szell I ^OOtT 
Mat subavashi. Makino fc Ebisuzaki 200 7;: [ Sesana et al.l 
[2008: 'Berentz en et al.Tl2009l : lAmaro-Seoane et al. I I2OI0I : 
[Scsana 2010) have been devoted to the study of MBHB 
dynamics in galactic nuclei, the major focus being the 
evolution of the binary semi-major axis to overcome the 
'last parsec problem'. Most of these works also report on 
the MBHB eccentricity evolution (which is in general more 
difficult to track and much more affected by numerical 
noise in A'^-body simulations), often ob serving a net in- 
crease during the hardening process (seei esana for 
a detailed discussion). This is of particular importance 
because (i) GW emission efficiency is a strong function 
of the eccentricity of the system (with eccentric binaries 
coalescing much faster) and (ii) even a small surviving 
eccentricity in the GW detection bands will be easily 
detectable l|Porter fc Sesana] I2OI0I ). possibly giving us 
clues about the binary evolution. Interestingly, the vast 
majority of the cited papers (with the notable exception of 
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iBerczik et"ani2006l : Ia maro- Seoane et al~ll20ld ) considered 
the MBHB evolution in non-rotating stellar systems. There 
are, however, three good reasons for considering rotating 
stellar distributions: 

(i) observationally, classical galaxy b ulges often show 
some degree of net rotation (see, e.g. iGadotti! l201ll l. 
and pseudobulges are mainly rotat ionally supported (e.g. 
iKormendv. Bender fc Cornell I [iooH ). A net rotation is also 
observed in the central region the Milky Way ijCenzel et al. I 
ll996l : ISchodel. Merritt fc Eckart Il2009( ). 

(ii) MBHBs formed during galaxy mergers are embed- 
ded in the remnant of the fusion of two galactic nuclei. Even 
assuming that the two nuclei originally had no spin, the 
orbital angular moment um associated with the merger will 
form a rotating system (iMilosavlievic fc Merritt1l200ll ). 

(iii) As noted bv llwasawa et aT counter- 
rotating stars are much more effective in extracting angu- 
lar momentum from the MBHB. We therefore expect the 
MBHB eccentricity to evolve differently depending on the 
degree of co(counter)-rotation of the surrounding stellar dis- 
tributions. 

In this letter, we study the eccentricity evolution of a 
MBHB with a small mass ratio (we consider q = M2/M1 — 
1/81) in stellar cusps with different degrees of rotation (from 
purely co-rotating to purely counter-rotating cusps). We ap- 
ply the hybrid formalism developed by Sesana et al. (2008, 
hereafter SHM08) that couples numerical three-body scat- 
terings to an analytical description of the cusp-binary inter- 
action. Given the several simplifying assumptions adopted in 
the hybrid model that may lead to spurious results when it 
comes to a delicate quantity like the binary eccentricity, we 
also performed calibrated direct A^'-body integrations of the 
binary-cusp system b y means of the di rect summation N- 
body code 0GRAPE jHarfst et al.ll2007l '). The letter is orga- 
nized as follows. In Section 2 we present a simple analytical 
argument to explain the different behaviours of co-rotating 
and counter-rotating stars interacting with the MBHBs, and 
we verify it with the support of calibrated three-body ex- 
periments. In Section 3 we describe the hybrid and TV-body 
models used to integrate the joint MBHB-cusp evolution. 
We present the results of the two methods, discussing sim- 
ilarities and differences in Section 4, and we conclude with 
some final remarks in Section 5. 



2 ANALYTICAL BACKGROUND AND THREE 
BODY SCATTERING 

In this section we discuss a simple argument to illustrate 
the different behaviour of co-rotating and counter-rotating 
stars in the star-binary interaction. For the sake of sim- 
plicity, we focus on a ideal coplanar case. A comprehen- 
sive analytical model will be presented in a follow up paper. 
We consider a system consisting of: (i) a binary with to- 
tal mass M = Mi + M2 (M2 <C Mi), semi-major axis a 
and eccentricity e, with initial energy £ = —GM\M2/{2a) 
and angular momentum C = Cz = ^^GMa{l — e^) (where 
fi = M-1M2/M) aligned along the positive z ajds; (ii) a star 
either co-rotating or counter-rotating with the binary , with 
m* ^ M2, in Keplerian orbit around Mi with semi- major 
axis a, « a and eccentricity e* ~ e. The star is charac- 
terized by an initial energy ~ ~GMim/{2a) and an- 



gular momentum C, — C^z ~ ±m«-\/GMia(l — e^) along 
the z axis (-1- if co-rotating with the binary, — if counter- 
rotating). Since M2 <S Mi, we ignore M2 in the energy 
and angular momentum budget of the stars, however such 
approximation (and the following dissertation) works fairly 
well also in the case of mildly unequal mass binaries with 
q — M2/M1 — 1/3. In any case, in the following we will con- 
sider Ml ~ M (thus, fi « M2 ). Setting a, « a and e, « e is 
particularly convenient for making a simple argument, since 
it cancels out complicated eccentricity dependencies. 

The starting point of our model is the definition of the 
MBHB eccentricity as a function of its energy £ and angular 
momentum magnitude £: 



2£C^ 
GM'^fiS 



Differentiation of equation ((T} leads to 

(l-e^/A£ 2A£. 
2e \ £ £z 
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2e 
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where x is defined by the last equality. We note that, being 
the binary angular momentum oriented along the z axis, 
and being in general A£ <^ jC, the eccentricity evolution 
depends on ACz only. Exchanges of the x and y component 
of jC will only result in a reorientation of the binary plane, 
without affecting e. The sign of Ae is therefore defined by 
the combination x of the energy and angular momentum 
variations. 

Following the three-body interaction, the star is ejected 
at infinity. The ejection is usually caused by a close en- 
counter with M2, that captures and ejects the star. This 
is known as the slingshot mechamsm. In general, counter- 
rotating stars have larger relative velocities with respect to 
AI2 then co-rotating ones, and therefore their capture and 
ejection cross sections are much smaller. In the presence of 
a binary, stars experience an asymmetric potential that al- 
lows for vari ations in the direction of their angular m omenta. 
As shown in iMerritt. Gualandris fc Mikkola I l|2009l ). an ec- 
centric binary exerts a semi-periodic forcing in a direction 
perpendicular to the orbital plane of the stars. Since the 
torque acting on the orbital plane of the stars is exerted 
by the binary, there is a correspondent change in the angu- 
lar momentum of the secondary, which in turn results in a 
change in the eccentricity. As a consequence of this torquing 
mechanism, both initially co-rotating and counter-rotating 
stars undergo secular evol ution, and are ejected when they 
co-rotate with the binary (jlwasawa et al. I2OI0I ). Following 
this observation, we assume the star's final energy to be neg- 
ligible (i.e. £,j, « 0) and its angular momentum z compo- 
nent to be positive, of the form C*,zf = rim,^GMia{l — e^) 
(where in general < < 1) . For co-rotating stars we there- 
fore have AjC^z = -ALz = m*{ri - l)^GMia{l - e^) and 



A£, 
tain 



-A£ = GMim,/(2a). Substituting in Eq. [2]we ob- 



Ae oc (2r; - 3) . 



(3) 



The same reasoning applies to the counter-rotating stars, 
with the exception that now the initial angular momentum 

s^), leading to 



is Ct 



,^GMia(l 

Ae oc {2ri + 1) 



(4) 
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Figure 1. Average star-binary exchanges in the three body inter- 
actions as a function of the MBHB eccentricity. Dashed curves: 
A£t^z = AZ2«_z/£J (i.e. the z component of the exchanged stel- 
lar angular momentum normalized to the angular momentum of 
a circular star with the same initial energy); dotted curves: * 
(i.e. normalized to the binding energy of a star with a* = a); 
solid curves x = x(-'^2 /'"■*)• Thick black curves are for counter- 
rotating stars, thin red curves are for co-rotating stars. 



This means that, unless ri > 3/2, the ejection of co-rotating 
stars decreases e and, conversely, the ejection of counter- 
rotating stars increases e. 

To test this simple heuristic description, we performed 
four sets of three-body scattering experiments. We consid- 
ered a binary with q = 1/81 and different eccentricities 
e = 0.1, 0.3, 0.6, 0.9. For each binary, we integrated the 
orbit of 2000 stars with a* ~ a, eccentricity drawn from a 
thermal distribution p(e) oc e (to mimic an isotropic distri- 
bution), and randomly oriented angular momentum £,,. We 
stored the energy and angular momentum of each star after 
the ejection, dividing the stars in two groups: the co-rotating 
(with > 0) and the counter-rotating (with £,,z < 0). 
Figure [T] shows the numerical results of the three-body ex- 
periments. The average star-binary energy exchange (dotted 
lines) is the same for both families and for any binary eccen- 
tricity. The angular momentum exchange is instead always 
much greater for counter-rotating stars, as expected. The 
practical consequence of this is that the quantity x oc Ae 
is positive for counter-rotating stars and negative for co- 
rotating stars at any given binary eccentricity. Our simple 
argument fixes e* — e, to cancel out complications due to 
eccentricity, and therefore can not catch the Ae dependence 
on e. However, the main result of our heuristic intuition is 
corroborated. 



3 INTEGRATION OF THE BINARY-STAR 
SYSTEM 

Having understood the different average behaviour of co- 
rotating and counter-rotating stars, we test here its prac- 
tical consequences. We follow the evolution of a MBHB in 
stellar cusps with different degrees of rotation. Integrations 
are performed both using an hybrid scheme and full direct 
A'^-body simulation. In the following, we briefly describe the 
two techniques. 



3.1 The hybrid model 

SHM08 constructed an hybrid model for evolving unequal 
MBHBs in stellar cusps. In short, it works as follows. They 
integrated 5 x 10* stars in bound orbits around A4i. In the 
Newtonian limit, the outcome of a scattering depends on 
a, /a = X, and can be scaled to any absolute value of a. 
Three-body scattering experiments provide the distributions 
of the average star-binary energy and angular momentum 
exchanges A£{x) AC{x), together with the bivariate distri- 
bution of the ejection times Af^x, t). Such information, ob- 
tained numerically, is then coupled to an analytical scheme 
for the joint evolution of the binary and the surrounding 
stellar distribution. For the practical integration of the hy- 
brid model, a stellar cusp with density p(r) oc r"'', normal- 
ized at the binary influence radius to an external isother- 
mal sphere (p(r) — cr^ / {2-KGr'^)), is assumed. The binary is 
placed with initial eccentricity at an initial separation 
where the mass in stars enclosed in its semi-major axis is 
twice the mass of M2 . In the hybrid model, stars at each x 
are weighted according to the radial density distribution, a 
detailed description of the technique is given in SHM08. 

The averaging procedure washes out the different be- 
havior of different types of stars having similar values of x. 
However, here we want to distinguish between co-rotating 
and counter-rotating stars. We therefore build two sets of 
distributions /S.£p/r(x), A.Cp/r(x), J\f^^^^{x,t), where p (pro- 
grade) denotes the average over co-rotating stars only, and 
r (retrograde) denotes the average over counter-rotating 
stars. Stellar cusps with different degrees of rotation are con- 
structed by mixing the two distribution as J-p -|- (1 — J-)r, 
where is a parameter running from zero to one, describ- 
ing the fraction of co-rotating stars. We have T — 0.5 for 
an isotropic cusp while fractions of 7^ > 0.5 (< 0.5) repre- 
sent cusps with a net degree of co-rotation (counter-rotation) 
with the binary. Note that, morphologically, the cusp is still 
spherically symmetric, i.e., this approximation does not take 
into account any possible axisymmetry or triaxiality induced 
in the cusp by the rotation. Moreover, being based on three- 
body scattering experiments, precession effects due to the 
extended potential of the cusp itself as well as other secular 
collective effects are ignored. 



3.2 N-body simulations 

In order to test the hybrid model, we also performed A^- 
body simulations of a MBHB embedded in a power-law 
stellar cusp and followed its evolution due to interactions 
with the stars. All the runs were performed wit h the direct- 
summation parallel A'^-body code <^GRAPE (' Harfst et al.l 
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Figure 2. Left: semi-major axis (upper panel) and eccentricity evolution (lower panel) of the MBHB according the the hybrid model. 
Different lines correspond to different fractions of co-rotating stars as labeled in the figure. We assume q = 1/81, 7 = —7/4 and = 0.5. 
Pi is the Keplerian orbital period of the binary at a^. Right: same for the direct Af-body integrations. Line styles as in the left plot. 



I2OO7I ) on the GPU enabled computers at the Max-Planck In- 
stitute for Astrophysics in Garching. The code uses a fourth- 
order Hermite scheme for the time-integration and can be 
used in combination with G RAPE or GPU hardw are, by 
means of the Sapporo library (|Gaburov et al.ir2009l ). 

While A^-body simulations are computationally much 
more expensive than three-body integrations and are there- 
fore limited in particle number, they allow to monitor the 
evolution of the black hole binary under the combined effects 
of star interactions, Newtonian precession and the Kozai 
mechanism. Precession of the orbits due to the distributed 
stellar mass is not taken into account in the hybrid model, 
but can play a role in the evolution of the binary. Moreover, 
the cumulative effect of the individual scatterings can sig- 
nificantly change the orientation of the bina ry orbital plane 
l|Merritt 11200^ : iGualandris fc Merritt 1120071 ). modifying the 
amount of cusp co-rotation as seen by the binary. 

The binary-stellar distribution model was constructed 
to match the initial set-up used in the hybrid scheme. 
We considered an unequal-mass binary with primary mass 
Ml = 10® Mq, mass ratio q — 1/81, initial semi-major axis 
Oi = 0.06 pc and initial eccentricity = 0.5. The stellar 
cusp follows a Bahcall-Wolf p(r) ~ |-~'^/^ profile at distances 
smaller than 1 pc, with total mass Mc ~ 2.5 X 10^ Mo and 
a mass enclosed in the binary orbit equal to 2M2. Because 
the distribution function adopted for the generation of the 
initial conditions is approximated, the constructed model is 
not in exact equilibrium. We therefore let the cusp relax 
before adding the secondary black hole. During this phase, 
the system undergoes a small expansion in the outer regions 
which however does not affect the distance range where the 
secondary black hole is placed. We use N = 32768 for all 
models which results in a black hole to star mass ratio of 
m./Mi = 7.5 X 10"^ 

In addition to a model with an isotropic distribution of 



velocities, we generate models with different fractions of co- 
rotating stars. These are obtained by reversing the sign of all 
the velocity components for a random subset of stars, at the 
time where the second black hole is added. This procedure is 
effectively equivalent to the mixing procedure Fp+{1 — T)r 
used in the hybrid scheme. 



4 RESULTS 

The main results of our experiments are collected in Fig. [2l 
where the MBHB semi-major axis and eccentricity evolution 
in cusps with different parameter is depicted. In all runs 
we used 7 = —7/4, q — 1/81 and — 0.5. The temporal 
evolution is plotted in units of the initial binary period Pi. 
The left plot in Fig.[2]shows evolutionary tracks produced by 
the hybrid integration scheme. The angular momentum of 
the star is non influential in the energy exchange, hence the 
orbital decay of the MBHB is hardly affected by the rotation 
of the cusp. The situation is drastically different in the case 
of the eccentricity evolution (lower panel). More counter- 
rotating cusps result in a faster evolution of the binary to- 
ward higher eccentricities. The critical fraction defining the 
transition between eccentricity growth and circularisation, 
for this particular choice of parameters, is ~ 0.7. Isotropic 
cusps (J^ = 0.5) lead to significant eccentricity growth, as 
found by SHM08. Results of direct A-body integrations are 
shown in the right panels of Fig. (2] for the same values of T 
used for the hybrid model. As expected, we find that cusps 
with a larger fraction of stars on co-rotating orbits with 
respect to the binary tend to circularise the binary while 
cusps with a larger fraction of stars on counter-rotating or- 
bits tend to increase the binary eccentricity, in agreement 
with the predictions from the hybrid model. The timescale 
for the binary decay is also in very good agreement with the 
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model's results. There are, however, some differences in the 
TV-body results with respect to the model. Firstly, the eccen- 
tricity in the A^-body integrations does not reach values as 
low and as high as in the hybrid model. In the fully counter- 
rotating cusp {J- = 0), the eccentricity grows to ~ 0.98 
while in the fully co-rotating cusp — 1) the eccentricity 
reaches ~ 0.078. This may be due to the small particle num- 
ber adopted in the simulations, but also to inaccuracies in 
the hybrid scheme when interpolations are performed to the 
eccentricity range boundaries. Also, the model with T — 0.5 
shows only a very slow growth in the eccentricity, in con- 
trast with what is found in the hybrid scheme. This may be 
the result of the suppression of the binary-induced secular 
evolution of the stars by Newtonian precession related to 
the extended cusp. As a consequence, counter-rotating stars 
are less prone to become co-rotating and are less efficient in 
extracting angular momentum from the MBHB. Lastly, the 
eccentricity growth in the counter-rotating T = case ap- 
pears much faster in the A^-body results than in the hybrid 
integration. This is due to the approximations used in the 
semi-analytic prescription. In the hybrid model, we subtract 
energy and angular momentum to the binary at the moment 
of the star ejection. Although this is a good approximation 
for the energy transfer, it is not for the angular momentum. 
Stars on counter-rotating orbits secularly subtract angular 
momentum to the binary to become co- rotating before being 
ejected. The MBHB eccentricity growth therefore does not 
occur at the moment of the star ejection (as in our hybrid 
scheme), but on a shorter timescale, during the star-MBHB 
interaction. This is reflected in the much faster eccentricity 
evolution in the counter-rotating N-body runs. 

As discussed in Section 3.2, the cumulative torque ex- 
erted by the stars changes also the orbital plane of the bi- 
nary. We find that the orientation changes by « 1 — 4 de- 
grees on the hardening timescale {a/a ~ lOOPo) in all the 
runs but the two most counter-rotating {T — 0,J- = 0.125). 
T his is in agreement wit h the change predicted in Eq. 4 
of iGualandris &: Merritt I l|2007l ) for an unequal mass bi- 
nary. We note that, in the counter-rotating models, the large 
change in the orbital plane (which results in the binary or- 
bital angular momentum reversal in the = case) occurs 
at t/Po > 200, i.e. after the bulk of the orbital evolution. 
This is because the large eccentricity attained by the binary 
(following the ejection of a large number of counter-rotating 
stars) results in a very small angular momentum. In this 
case, small torques can easily produce a drastic reorienta- 
tion of the binary plane. This dynamical aspect does not 
affect our results, but is worth further investigation. 



5 FINAL REMARKS 

We demonstrate that the degree of rotation of the stel- 
lar background surrounding a MBHB determines the ec- 
centricity evolution of the binary. This has already been 
discussed for wide BH pairs, subje ct to dynamical fr iction 
exerted by a rotating background l|Dotti et al. Il2007t) . and 
in the e arly evolution of pairing M BHBs in rotating star 
clusters (|Amaro-Seoane et al. II2OI0I ). This latter work, in 
particular, describes the different action of dynamical fric- 
tion in counter-rotating clusters, and it is somewhat comple- 
mentary to our paper. Here we highlight for the first time 



the physical mechanism at work in close binaries evolving 
through interactions with single stars. We find that, for an 
unequal mass binary, stellar systems co-rotating with re- 
spect to the binary tend to circularise its orbit. On the 
other hand, if stars have, on average, angular momenta anti- 
aligned with respect to the binary, a strong increase in the 
eccentricity is observed. 

The dependence of the eccentricity evolution of the bi- 
nary on the rotation of the stellar background received so far 
little attention. Rotation can be due to secular evolution or 
merger events. For equal mass binaries, presumably formed 
through major galaxy mergers, the MBHs are likely to co- 
rotate with the nucleus of the remnant, reminiscent of the 
orbital motion of the parent nuclei. In this case, more circu- 
lar MBHBs are e xpected. Very unequal mas s binaries could 
form in situ (see I Amaro-Seoane et al." 1 l2007l . and references 
therein), or via galactic minor mergers. In this case, the 
original rotation of the nucleus of the primary galaxy would 
be less perturbed by the interaction, possibly resulting in 
counter-rotating systems, and extremely eccentric binaries. 
In both equal and unequal cases, the exact eccentricity evo- 
lution depends on the degree of rotation and on the slope of 
the stellar profile. A more detailed analysis is postponed to 
a future investigation. 
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